home *** CD-ROM | disk | FTP | other *** search
/ Hottest 6 / Hottest 6 (1996)(PDSoft)[!].iso / software / videoutils / a-g / ami2d / test / srm4.rexx < prev    next >
OS/2 REXX Batch file  |  1978-11-24  |  11KB  |  481 lines

  1. /* simple solid rocket motor */
  2. options results
  3.  
  4. if ~show('ports', "AMI2D") then do
  5.     address command 'run ami2d:bin/ami2d'
  6.     do while ~show('ports', "AMI2D")
  7.     end
  8. end
  9.  
  10. address ami2d
  11.  
  12. call 'ami2d:rexx/pstrain'
  13.  
  14. if ~show('l', "rexxmathlib.library") then do
  15.     check = addlib('rexxmathlib.library',0,-30,0)
  16. end
  17.  
  18. pi = 3.141592654
  19. rb = 6.0
  20. ra = 1.0
  21. rs = 4.8
  22. w2 = 0.2
  23. h = 0.05
  24. E = 1500.0
  25. nu = 0.4900
  26. alp = 6.0e-6
  27. Ec = 10.0e6
  28. nuc = 0.3000
  29. alpc = 5.0e-5
  30. P = -1000.0
  31. Tref = 140
  32. To = 77
  33. T = 500
  34.  
  35. 'reset'
  36. 'problem(thermal)'
  37. 'iso(1,'E','nu',thermal,'alp','Tref')'
  38. 'iso(2,'Ec','nuc',thermal,'alpc','Tref')'
  39. 'therm(1,.5)'
  40. 'therm(2,1)'
  41. 'coord(cart)'
  42. a1 = asin(w2/ra)
  43. a2 = pi/4 - a1
  44. da = (a2 - a1)/16
  45. b1 = asin(w2/rs)
  46. b2 = pi/4 - b1
  47. db = (b2 - b1)/16
  48. x1 = ra*cos(a1)
  49. y1 = ra*sin(a1)
  50. x2 = rs
  51. y2 = w2
  52. g = 1.5
  53. i = line(x1,y1,x2,y2,g,1,21,1)
  54. xc = rs
  55. yc = 0.0
  56. i = arc(xc,yc,pi/2,0.0,w2,21,9,1)
  57. x1 = ra*cos(a1 + da)
  58. y1 = ra*sin(a1 + da)
  59. x2 = rs*cos(b1 + db)
  60. y2 = rs*sin(b1 + db)
  61. i = line(x1,y1,x2,y2,g,41,11,2)
  62. i = ltip(xc,yc,x2,y2,0.0,0.25,w2,rb,61,3,2)
  63. x1 = ra*cos(a1 + 2*da)
  64. y1 = ra*sin(a1 + 2*da)
  65. x2 = rs*cos(b1 + 2*db)
  66. y2 = rs*sin(b1 + 2*db)
  67. i = line(x1,y1,x2,y2,g,81,21,1)
  68. i = ltip(xc,yc,x2,y2,0.0,0.5,w2,rb,101,5,1)
  69. x1 = ra*cos(a1 + 3*da)
  70. y1 = ra*sin(a1 + 3*da)
  71. x2 = rs*cos(b1 + 3*db)
  72. y2 = rs*sin(b1 + 3*db)
  73. i = line(x1,y1,x2,y2,g,121,11,2)
  74. i = ltip(xc,yc,x2,y2,0.0,0.75,w2,rb,141,3,2)
  75. x1 = ra*cos(a1 + 4*da)
  76. y1 = ra*sin(a1 + 4*da)
  77. x2 = rs*cos(b1 + 4*db)
  78. y2 = rs*sin(b1 + 4*db)
  79. i = line(x1,y1,x2,y2,g,161,21,1)
  80. i = ltip(xc,yc,x2,y2,0.0,1.0,w2,rb,181,5,1)
  81. i = arc(0.0,0.0,pi/16,0.0,rb + h/2,190,3,2)
  82. i = arc(0.0,0.0,pi/16,0.0,rb + h,195,5,1)
  83.  
  84. x1 = ra*cos(a1 + 6*da)
  85. y1 = ra*sin(a1 + 6*da)
  86. x2 = rs*cos(3*pi/32)
  87. y2 = rs*sin(3*pi/32)
  88. i = line(x1,y1,x2,y2,g,201,11,2)
  89. x1 = x2
  90. y1 = y2
  91. x2 = rb*cos(3*pi/32)
  92. y2 = rb*sin(3*pi/32)
  93. i = line(x1,y1,x2,y2,1.0,221,3,2)
  94. x1 = x2
  95. y1 = y2
  96. x2 = (rb + h)*cos(3*pi/32)
  97. y2 = (rb + h)*sin(3*pi/32)
  98. i = line(x1,y1,x2,y2,1.0,225,2,2)
  99.  
  100. x1 = ra*cos(a1 + 8*da)
  101. y1 = ra*sin(a1 + 8*da)
  102. x2 = rs*cos(pi/8)
  103. y2 = rs*sin(pi/8)
  104. i = line(x1,y1,x2,y2,g,241,21,1)
  105. x1 = x2
  106. y1 = y2
  107. x2 = rb*cos(pi/8)
  108. y2 = rb*sin(pi/8)
  109. i = line(x1,y1,x2,y2,1.0,261,5,1)
  110. x1 = x2
  111. y1 = y2
  112. x2 = (rb + h)*cos(pi/8)
  113. y2 = (rb + h)*sin(pi/8)
  114. i = line(x1,y1,x2,y2,1.0,265,3,1)
  115. x1 = ra*cos(a1 + 10*da)
  116. y1 = ra*sin(a1 + 10*da)
  117. x2 = rs*cos(5*pi/32)
  118. y2 = rs*sin(5*pi/32)
  119. i = line(x1,y1,x2,y2,g,281,11,2)
  120. x1 = x2
  121. y1 = y2
  122. x2 = rb*cos(5*pi/32)
  123. y2 = rb*sin(5*pi/32)
  124. i = line(x1,y1,x2,y2,1.0,301,3,2)
  125. x1 = x2
  126. y1 = y2
  127. x2 = (rb + h)*cos(5*pi/32)
  128. y2 = (rb + h)*sin(5*pi/32)
  129. i = line(x1,y1,x2,y2,1.0,305,2,2)
  130.  
  131. xc = rs*cos(pi/4)
  132. yc = rs*sin(pi/4)
  133. x1 = ra*cos(a1 + 12*da)
  134. y1 = ra*sin(a1 + 12*da)
  135. x2 = rs*cos(b1 + 12*db)
  136. y2 = rs*sin(b1 + 12*db)
  137. i = line(x1,y1,x2,y2,g,321,21,1)
  138. i = rtip(xc,yc,x2,y2,pi/4,1.0,w2,rb,341,5,1)
  139. i = arc(0.0,0.0,3*pi/16,pi/4,rb + h/2,350,3,2)
  140. i = arc(0.0,0.0,3*pi/16,pi/4,rb + h,355,5,1)
  141.  
  142. x1 = ra*cos(a1 + 13*da)
  143. y1 = ra*sin(a1 + 13*da)
  144. x2 = rs*cos(b1 + 13*db)
  145. y2 = rs*sin(b1 + 13*db)
  146. i = line(x1,y1,x2,y2,g,361,11,2)
  147. i = rtip(xc,yc,x2,y2,pi/4,0.75,w2,rb,381,3,2)
  148.  
  149. x1 = ra*cos(a1 + 14*da)
  150. y1 = ra*sin(a1 + 14*da)
  151. x2 = rs*cos(b1 + 14*db)
  152. y2 = rs*sin(b1 + 14*db)
  153. i = line(x1,y1,x2,y2,g,401,21,1)
  154. i = rtip(xc,yc,x2,y2,pi/4,0.5,w2,rb,421,5,1)
  155.  
  156. x1 = ra*cos(a1 + 15*da)
  157. y1 = ra*sin(a1 + 15*da)
  158. x2 = rs*cos(b1 + 15*db)
  159. y2 = rs*sin(b1 + 15*db)
  160. i = line(x1,y1,x2,y2,g,441,11,2)
  161. i = rtip(xc,yc,x2,y2,pi/4,0.25,w2,rb,461,3,2)
  162.  
  163. x1 = ra*cos(a2)
  164. y1 = ra*sin(a2)
  165. x2 = xc + w2*sin(pi/4)
  166. y2 = yc - w2*cos(pi/4)
  167. g = 1.5
  168. i = line(x1,y1,x2,y2,g,481,21,1)
  169. i = arc(xc,yc,-pi/4,pi/4,w2,501,9,1)
  170.  
  171. a1 = asin(w2/ra) + pi/4
  172. a2 = pi/2 - asin(w2/ra)
  173. da = (a2 - a1)/16
  174. b1 = asin(w2/rs) + pi/4
  175. b2 = pi/2 - asin(w2/rs)
  176. db = (b2 - b1)/16
  177.  
  178. i = arc(xc,yc,3*pi/4,pi/4,w2,541,9,1)
  179. x1 = ra*cos(a1)
  180. y1 = ra*sin(a1)
  181. x2 = xc - w2*sin(pi/4)
  182. y2 = yc + w2*cos(pi/4)
  183. g = 1.5
  184. i = line(x1,y1,x2,y2,g,521,21,1)
  185. x2 = rs*cos(b1 + db)
  186. y2 = rs*sin(b1 + db)
  187. i = ltip(xc,yc,x2,y2,pi/4,0.25,w2,rb,581,3,2)
  188. x1 = ra*cos(a1 + da)
  189. y1 = ra*sin(a1 + da)
  190. i = line(x1,y1,x2,y2,g,561,11,2)
  191. x2 = rs*cos(b1 + 2*db)
  192. y2 = rs*sin(b1 + 2*db)
  193. i = ltip(xc,yc,x2,y2,pi/4,0.5,w2,rb,621,5,1)
  194. x1 = ra*cos(a1 + 2*da)
  195. y1 = ra*sin(a1 + 2*da)
  196. i = line(x1,y1,x2,y2,g,601,21,1)
  197. x2 = rs*cos(b1 + 3*db)
  198. y2 = rs*sin(b1 + 3*db)
  199. i = ltip(xc,yc,x2,y2,pi/4,0.75,w2,rb,661,3,2)
  200. x1 = ra*cos(a1 + 3*da)
  201. y1 = ra*sin(a1 + 3*da)
  202. i = line(x1,y1,x2,y2,g,641,11,2)
  203. x2 = rs*cos(b1 + 4*db)
  204. y2 = rs*sin(b1 + 4*db)
  205. i = ltip(xc,yc,x2,y2,pi/4,1.0,w2,rb,701,5,1)
  206. i = arc(0.0,0.0,5*pi/16,pi/4,rb + h/2,710,3,2)
  207. i = arc(0.0,0.0,5*pi/16,pi/4,rb + h,715,5,1)
  208. x1 = ra*cos(a1 + 4*da)
  209. y1 = ra*sin(a1 + 4*da)
  210. i = line(x1,y1,x2,y2,g,681,21,1)
  211.  
  212. x1 = rs*cos(11*pi/32)
  213. y1 = rs*sin(11*pi/32)
  214. x2 = rb*cos(11*pi/32)
  215. y2 = rb*sin(11*pi/32)
  216. i = line(x1,y1,x2,y2,1.0,741,3,2)
  217. x1 = x2
  218. y1 = y2
  219. x2 = (rb + h)*cos(11*pi/32)
  220. y2 = (rb + h)*sin(11*pi/32)
  221. i = line(x1,y1,x2,y2,1.0,745,2,2)
  222. x1 = ra*cos(a1 + 6*da)
  223. y1 = ra*sin(a1 + 6*da)
  224. x2 = rs*cos(11*pi/32)
  225. y2 = rs*sin(11*pi/32)
  226. i = line(x1,y1,x2,y2,g,721,11,2)
  227.  
  228. x1 = rs*cos(3*pi/8)
  229. y1 = rs*sin(3*pi/8)
  230. x2 = rb*cos(3*pi/8)
  231. y2 = rb*sin(3*pi/8)
  232. i = line(x1,y1,x2,y2,1.0,781,5,1)
  233. x1 = x2
  234. y1 = y2
  235. x2 = (rb + h)*cos(3*pi/8)
  236. y2 = (rb + h)*sin(3*pi/8)
  237. i = line(x1,y1,x2,y2,1.0,785,3,1)
  238. x1 = ra*cos(a1 + 8*da)
  239. y1 = ra*sin(a1 + 8*da)
  240. x2 = rs*cos(3*pi/8)
  241. y2 = rs*sin(3*pi/8)
  242. i = line(x1,y1,x2,y2,g,761,21,1)
  243.  
  244. x1 = rs*cos(13*pi/32)
  245. y1 = rs*sin(13*pi/32)
  246. x2 = rb*cos(13*pi/32)
  247. y2 = rb*sin(13*pi/32)
  248. i = line(x1,y1,x2,y2,1.0,821,3,2)
  249. x1 = x2
  250. y1 = y2
  251. x2 = (rb + h)*cos(13*pi/32)
  252. y2 = (rb + h)*sin(13*pi/32)
  253. i = line(x1,y1,x2,y2,1.0,825,2,2)
  254. x1 = ra*cos(a1 + 10*da)
  255. y1 = ra*sin(a1 + 10*da)
  256. x2 = rs*cos(13*pi/32)
  257. y2 = rs*sin(13*pi/32)
  258. i = line(x1,y1,x2,y2,g,801,11,2)
  259.  
  260. xc = 0.0
  261. yc = rs
  262. x2 = rs*cos(b1 + 12*db)
  263. y2 = rs*sin(b1 + 12*db)
  264. i = rtip(xc,yc,x2,y2,pi/2,1.0,w2,rb,861,5,1)
  265. i = arc(0.0,0.0,7*pi/16,pi/2,rb + h/2,870,3,2)
  266. i = arc(0.0,0.0,7*pi/16,pi/2,rb + h,875,5,1)
  267. x1 = ra*cos(a1 + 12*da)
  268. y1 = ra*sin(a1 + 12*da)
  269. i = line(x1,y1,x2,y2,g,841,21,1)
  270.  
  271. x2 = rs*cos(b1 + 13*db)
  272. y2 = rs*sin(b1 + 13*db)
  273. i = rtip(xc,yc,x2,y2,pi/2,0.75,w2,rb,901,3,2)
  274. x1 = ra*cos(a1 + 13*da)
  275. y1 = ra*sin(a1 + 13*da)
  276. i = line(x1,y1,x2,y2,g,881,11,2)
  277. x2 = rs*cos(b1 + 14*db)
  278. y2 = rs*sin(b1 + 14*db)
  279. i = rtip(xc,yc,x2,y2,pi/2,0.5,w2,rb,941,5,1)
  280. x1 = ra*cos(a1 + 14*da)
  281. y1 = ra*sin(a1 + 14*da)
  282. i = line(x1,y1,x2,y2,g,921,21,1)
  283. x2 = rs*cos(b1 + 15*db)
  284. y2 = rs*sin(b1 + 15*db)
  285. i = rtip(xc,yc,x2,y2,pi/2,0.25,w2,rb,981,3,2)
  286. x1 = ra*cos(a1 + 15*da)
  287. y1 = ra*sin(a1 + 15*da)
  288. i = line(x1,y1,x2,y2,g,961,11,2)
  289.  
  290. i = arc(xc,yc,0.0,pi/2,w2,1021,9,1)
  291. x1 = ra*cos(a2)
  292. y1 = ra*sin(a2)
  293. x2 = w2
  294. y2 = rs
  295. g = 1.5
  296. i = line(x1,y1,x2,y2,g,1001,21,1)
  297.  
  298. do j=1 to 14
  299.     n = j*2 - 1
  300.     'quad8('j',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  301.     'press('j',1,'P')'
  302. end
  303. do j=1 to 14
  304.     n = j*2 + 79
  305.     'quad8('j + 20',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  306. end
  307. 'quad8(35,2,187,189,199,197,188,194,198,192)'
  308. 'quad8(36,2,185,187,197,195,186,192,196,190)'
  309. 'group(case,add,185,199)'
  310. do j=1 to 12
  311.     n = j*2 + 159
  312.     'quad8('j + 40',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  313. end
  314. 'quad8(53,2,265,185,195,267,225,190,227,266)'
  315. 'group(case,add,225,227)'
  316. 'group(case,add,265,267)'
  317. do j=1 to 12
  318.     n = j*2 + 239
  319.     'quad8('j + 60',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  320. end
  321. 'quad8(73,2,345,265,267,355,305,266,307,350)'
  322. 'group(case,add,305,307)'
  323. do j=1 to 14
  324.     n = j*2 + 319
  325.     'quad8('j + 80',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  326. end
  327. 'quad8(95,2,347,345,355,357,346,350,356,352)'
  328. 'quad8(96,2,349,347,357,359,348,352,358,354)'
  329. 'group(case,add,345,359)'
  330. do j=1 to 14
  331.     n = j*2 + 399
  332.     'quad8('j + 100',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  333.     'press('j + 100',3,'P')'
  334. end
  335.  
  336. do j=1 to 13
  337.     n = j*2 + 519
  338.     'quad8('j + 120',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  339.     'press('j + 120',1,'P')'
  340. end
  341. 'quad8(134,1,547,509,429,627,548,469,628,587)'
  342. 'press(134,1,'P')'
  343. do j=1 to 13
  344.     n = j*2 + 599
  345.     'quad8('j + 140',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  346. end
  347. 'quad8(154,1,627,429,349,707,628,389,708,667)'
  348. 'quad8(155,2,707,349,359,717,708,354,718,712)'
  349. 'quad8(156,2,705,707,717,715,706,712,716,710)'
  350. 'group(case,add,705,718)'
  351. do j=1 to 12
  352.     n = j*2 + 679
  353.     'quad8('j + 160',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  354. end
  355. 'quad8(173,2,785,705,715,787,745,710,747,786)'
  356. 'group(case,add,745,747)'
  357. 'group(case,add,785,787)'
  358. do j=1 to 12
  359.     n = j*2 + 759
  360.     'quad8('j + 180',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  361. end
  362. 'quad8(193,2,865,785,787,875,825,786,827,870)'
  363. 'group(case,add,825,827)'
  364. do j=1 to 14
  365.     n = j*2 + 839
  366.     'quad8('j + 200',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  367. end
  368. 'quad8(215,2,867,865,875,877,866,870,876,872)'
  369. 'quad8(216,2,869,867,877,879,868,872,878,874)'
  370. 'group(case,add,865,879)'
  371. do j=1 to 14
  372.     n = j*2 + 919
  373.     'quad8('j + 220',1,'n','n + 2','n + 82','n + 80','n + 1','n + 42','n + 81','n + 40')'
  374.     'press('j + 220',3,'P')'
  375. end
  376. do j=1 to 221 by 20
  377.     'press('j',4,'P')'
  378. end
  379. do j=29 to 189 by 40
  380.     'dispy('j',0.0)'
  381. end
  382. 'dispy(194,0.0)'
  383. 'dispy(199,0.0)'
  384. do j=869 to 1029 by 40
  385.     'dispx('j',0.0)'
  386. end
  387. 'dispx(874,0.0)'
  388. 'dispx(879,0.0)'
  389. do j=1 to 29
  390.     'temp('j','T')'
  391. end
  392. do j=481 to 509
  393.     'temp('j','T')'
  394. end
  395. do j=521 to 548
  396.     'temp('j','T')'
  397. end
  398. do j=1001 to 1029
  399.     'temp('j','T')'
  400. end
  401. do j=41 to 441 by 40
  402.     'temp('j','T')'
  403. end
  404. do j=561 to 961 by 40
  405.     'temp('j','T')'
  406. end
  407. do j=195 to 199
  408.     'temp('j','To')'
  409. end
  410. 'temp(227,'To')'
  411. 'temp(267,'To')'
  412. 'temp(307,'To')'
  413. do j=355 to 359
  414.     'temp('j','To')'
  415. end
  416. do j=715 to 718
  417.     'temp('j','To')'
  418. end
  419. 'temp(747,'To')'
  420. 'temp(787,'To')'
  421. 'temp(827,'To')'
  422. do j=875 to 879
  423.     'temp('j','To')'
  424. end
  425. exit
  426.  
  427. line:
  428. arg xa,ya,xb,yb,gg,i1,n,di
  429.  
  430. n1 = di*(n - 1)
  431. do i=0 to n1  by di
  432.     j = i1 + i
  433.     xn = xa + (xb-xa)*pow(i/n1,gg)
  434.     yn = ya + (yb-ya)*pow(i/n1,gg)
  435.     'node('j','xn','yn')'
  436. end
  437. return j
  438.  
  439. arc:
  440. arg cx,cy,aa,ab,r1,i1,n,di
  441.  
  442. n1 = di*(n - 1)
  443. do i=0 to n1  by di
  444.     j = i1 + i
  445.     an = aa + (ab-aa)*i/n1
  446.     xn = r1*cos(an) + cx
  447.     yn = r1*sin(an) + cy
  448.     'node('j','xn','yn')'
  449. end
  450. return j
  451.  
  452. ltip:
  453. arg xc,yc,x1,y1,ab,fn,w2,rb,i1,n,di
  454.  
  455. xa = w2*cos(ab + pi/4) + xc
  456. ya = w2*sin(ab + pi/4) + yc
  457. xb = rb*cos(ab + pi/16)
  458. yb = rb*sin(ab + pi/16)
  459. xx = xa + (xb-xa)*fn
  460. yy = ya + (yb-ya)*fn
  461. i = line(x1,y1,xx,yy,1.0,i1,n,di)
  462. rr = sqrt(xx*xx + yy*yy)
  463. aa = atan(yy/xx)
  464. i = arc(0.0,0.0,aa,ab,rr,i,n,di)
  465. return i
  466.  
  467. rtip:
  468. arg xc,yc,x1,y1,ab,fn,w2,rb,i1,n,di
  469.  
  470. xa = w2*cos(ab - pi/4) + xc
  471. ya = w2*sin(ab - pi/4) + yc
  472. xb = rb*cos(ab - pi/16)
  473. yb = rb*sin(ab - pi/16)
  474. xx = xa + (xb-xa)*fn
  475. yy = ya + (yb-ya)*fn
  476. i = line(x1,y1,xx,yy,1.0,i1,n,di)
  477. rr = sqrt(xx*xx + yy*yy)
  478. aa = atan(yy/xx)
  479. i = arc(0.0,0.0,aa,ab,rr,i,n,di)
  480. return i
  481.